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Abstract 

In this paper we present the discontinuous Galerkin method to solve the 
problem of the two-dimensional air pollution model. The resulting system 
of ordinary differential equations is called the semidiscrete formulation. We 
show the existence and uniqueness of the ODE system and provide the error 
estimates for the numerical error. 
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1. Introduction 

Air pollution is the introduction of chemicals, particulate matter, or bi- 
ological materials that cause harm or discomfort to humans or other living 
organisms, or cause damage to the natural environment or built environ- 
ment, into the atmosphere. The basic technology for analyzing air pollution 
is through the mathematical models and numerical methods for predicting 
the transport of air pollutants in the lower atmosphere [H, 0, 0, Isj . Different 
air pollution models have been developed in the last decades by the National 
Environmental Research Institute (http : //www . dmu . dk/ en/air/models/) . 
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In the present paper we consider the following Danish Eulerian model j2|, 0, Isf 



f{u) = -{h + k2)u + E + Q{u), (lb) 
?/, 0) = uoix, y), {x, y) G fi, (Ic) 
t) Ian = 0, t G [0, T]. (Id) 

The different quantities involved in the mathematical model have the 
following meaning: 

• the concentration is denoted by u] 

• c and e are wind velocities; 

• kx and ky are diffusion coefficients; 

• the emission source is described by E; 

• ki and k2 are constant deposition coefficients; 

• the chemical reaction is denoted by Q. 
Meanwhile, we give the following assumptions: 

• ue Hlin) n H^in), ut, uu e L'^in)- 

• Q{u) satisfy the Lipschitz condition; 

• < fc* ^ min{|/c2;| , \ky\} ^ max{|/c2;| , \ky\} ^ k*, < c* ^ 
min{|c| , |e|} ^ max{|c| , |e|} ^ c*, fc*, k*, c*, c*are constants. 

A general description of the Danish Eulerian Model and its numerical 
treatment is given in j^, @, 0]- Research on the finite difference method and 
finite volume element method for this air pollution model already has good 
results [^,@, l3, 11 1 . In this article, we use the discontinuous Galerkin method 



(DG method) to analyse and solve the air pollution model. 

DG methods in mathematics form a class of numerical methods for solving 
partial differential equations. They have recently gained popularity due to 
many of their attractive properties, refer to 12, 13, 14, 15, 16, 13, E, [isl, 20 
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First of all, the flexibility of the methods allows for general non-conforming 
meshes with variable degree of approximation. This makes the implementa- 
tion of h-p adaptivity for DG easier than that for conventional approaches. 
Moreover, the DG methods are locally mass conservative at the element level. 
In addition, they have less numerical diffusion than most conventional algo- 
rithms, thus are likely to offer more accurate solution for at least advection- 
dominated transport problems. They handle rough coefficient problems and 
capture the discontinuity in the solution very well by the nature of discontin- 
uous function space. Furthermore, the DG methods are easier to implement 
than most traditional finite element methods. The trial and test spaces are 
easier to construct than conforming methods because they are local. 

The paper is organized as follows: In Section 2, the variational formula- 
tion of the DG method is stated. And we show the existence and uniqueness 
of the resulting ordinary differential equations system. Finally we provide 
the error estimates for the numerical error in Section 3. 

2. Semidiscrete formulation 

In this section, we approximate the solution u{t) by a function Uh{t) that 
belongs to the finite-dimensional space I)k{eh) for all t > 0. The solution 
Uh is referred to as the semidiscrete solution. In what follows, we assume 
that s > |. We introduce a bilinear form Jg"'^": H^^Sh) x H^{eh) M that 
penalize the jump of the function values: 



The parameter is called penalty parameter. It is nonnegative real number. 
The power /3o is positive number. |e| simply means the length of e. We now 
define the DG bihnear form : H^{eh) x H^{eh) — )■ M 



The bilinear form contains another parameter e that may take the value 
-1,0, or 1. is symmetric if e = — 1 and it is nonsymmetric otherwise. 




a, 
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This bilinear form yields the following energy seminorm: 



Second, the convection term ■§^{cu) + ^(eu) is approximated by an upwind 
discretization. Let us denote the upwind value of a function w by w"^. We 

recall that ^ ^ ^ is a unit normal vector pointing from to El: 

{If I ^1 , if cnt + ent > 10 
■=^4 4~ ^e = dElndEl. 
w\ei, if crii + en2 < 

Let 

6(c, e; w, t;) = - ^ / {cw— + a") + / (cr^iW^^'H + en2w"*'H) 

The general semidiscrete DG variational formulation of problem ( |Ta|) -( ]ld[l 
IS as follows: Find Uh G L^{0,T;Vk{eh)), such that 

dU 

Wt >0,yve Vkieh), i-^, v)n + a,{Uh{t),v) + 6(c, e; [/^(t), t;) = L(t/;,(t), t;), 

(2a) 

Vt; e (f//,(0), = (mo, t;)f,, (2b) 

where the form L is 

L(w;t;) = / f{w)v. 



The next lemma establishes the consistency between the model problem 
and the variational formulation. 



Lemma 2.1. Assume that the weak solution u of problem pm) - pd) belongs 
to H^{0,T; H'^{6h)) , then u satisfies the variational problem [2c\) - ^2B\] . 

Proof. Let f be a test function in Vk{eh). We multiply by t>|£; and integrate 
by parts on one element E G 6h, and use Green's theorem: 

.du r ^ du dv ^ dudv dv dv 

dt^ ^ Je ^ dx dx ^ dy dy dx dy 



f , , du^ , du^ V _v , f 

/ (—Kj.-—niV — kj,—n2V + cuniV + eun2V)= i 
JdE dx ''dy Je 



f{u)v 
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Summing over all elements and using the regularity of the exact solution, we 
obtain 

, du r ^ dudv ^ dudv dv dv 

dt ' ^ ^ Je ^ dx dx ^ dy dy dx dy 

/^^(^ ^^^^ 
{-{kx-Q^ni}[v] - {ky-g^nt}[v] + cuni[v] + eun^[v]) + 



f{u)v. 



Since Uup = u, we clearly have our result. □ 

2.1. Existence and uniqueness of the solution 

Because of the lack of continuity constraints between mesh elements for 
the test functions, the basic functions of 'Dk{sh)) have a support contained 
in one element. We write 

'Dk{eh) = span{(t)f : 1 ^ i ^ Nioc, E e Eh} 
with _ 

,£/ X \ (1^1° Fe{x), X e e, 

0, X f E. 

In 2D, we have (p{x,y) — x^y^,I + J = i,0 < i < k. This yields the local 
dimension 

{k+l){k + 2) 
^^loc — 2 ■ 

using the global basis functions, we can expand the semidiscrete solution 

Nloc 

Vt e (0, T), y{x, y) e Unit, x,y)=Y, Y.^fit)^f{^^ v)- (3) 

The degrees of freedom ^^^s are functions of time. Let Nf.i denote the number 
of elements in the mesh. We can rename the basis functions and the degrees 
of freedom such that 

{0f : 1 < i < Nloc, Eeeh} = {4>j : 1 < J < iV,„,iVeJ, 
{if:l<i< Nloc, Eesh} = {Cj : 1 < J < NiocNd}. 
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Plugging ([3]) into the variational problem fl2ap - pb|) yields a linear system of 
ordinary differential equations as follows: 

M^(t) + (A + i?)^=G(0, 
Mm = Uo. 

The matrices M, A are called the mass and stiffness matrices, and they are 
defined by 

VI ^i, j ^ NiocNei, Mij = {(t>j, (j)i)n, Aij = a^{(f)j, (pi). 

The matrix B results from the convective term, and the vector G{^) 
depends on the vector solution 

VI ^i, j ^ NiocNei, {B)ij = b{c; e; (f)j, (f)i), 

Since the matrix M is invertible and the vector function G{S,) is Lipschitz 
with respect to ^, there exists a unique solution to the variational problem 
(l2iD-(l2b]). 



3. Error estimates 

In this section, we first present the Gronwall's inequalities [2l|, which are 
important tools for analyzing time-dependent problems. 

Lemma 3.1 (Continuous Gronwall inquality). Let /, (7, h he piecewise contin- 
uous nonnegative functions defined on (a, b). Assume that g is nondereasing . 
Assume that there is a positive constant C independent of t such that 

Vt G (a, b), fit) + hit) < git) + C I fis)ds. 

J a 

Then, 

Vt G (a, 6), fit) + hit) < e^(*-")^(t). 



Now we state a priori error estimates for the semidiscrete scheme [22 
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Theorem 3.2. Assume that the solution u to problem ( f7a|) -/ l7^) belongs to 
H^{0,T; H'^{eh)) and that Uq belongs to H'^{eh) for s > 3/2. Assume that 
/3o > 1. In the case of SIPG and IIPG, assume that is sufficiently large 
for all e. Then, there is a constant C independent of h such that 

f r 2 ^ 

+ / \\u{t)-UHm;dt 







(eh)' 



Proof. We omit some details which is similar to the proof of Theorem 3. 4. ([13!]) 
We write u — Uh = p — x with p = u — u and x = Uh — u. The function 
u G 'Dk{eh)) is an approximation of u that satisfies good error bounds. The 
error equation is satisfied for all v in Vk{eh)- 

dx dp 
(-^, v)n + a^ix, v) + b{c, e; x, v) = (-^, v)n + a^{p, v) 

+ b{c,e;p,v) + {fiUH)-f{u),v)n. 
Now, by choosing v = x and using the coercivity of a^, we obtain 

^ ^ II l|2 , II ||2 , u \ ^ (^P \ 

2Jt WxWmn) + 1^ Me + b{c; e; X, X) ^ X)n 
+ aeip, x) + &(c, e; p, x) + ifiUh) - /(«), x)n- 

We use Green's formula and the fact that V ■ ( ^ 1=0: 



h 
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Thus we obtain 



5^ /(c^ + e^)(x"1x]-i[x1) 
5^ /'(c^ + e^)(x"^[x]-{x}M) 



We now bound each term in b{c, e; p, x)- Using Cauchy-Schwarz's and Young's 
inequahties, we have 



l|VxlL.(i,)< g 11x11. + 

E&Sh 



and 



1 

4 



0, e 
2 



\cni + e?22rP* 



0, e 



0, e 



r* II n"P|l^ 

^ IIP llL2(e). 

Finally, we bound the nonlinear source term, using the Lipschitz property: 

ifiU,) - f{u))x < C\\{U, - ^^)\\LHn)\\xLHn) < CMl^^n) + CMlnny 

The other terms are identical to the ones in the proof of Theorem 2.13 and 
3.4 ([13!]). Then the main result is obtained by combining all bounds and 
using Gronwall's inequality of Lemma 13.11 □ 

We can choose any of the time discretizations such as backward Euler and 
forward Euler and some that are of high order such as Crank-Nicolson and 
Runge-Kutta methods. The analysis of the resulting fully discrete schemes 
can be done in a common way. 
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